knitr::opts_chunk$set(echo = TRUE)

cmdstanr::set_cmdstan_path(path = "C:/Users/kueng/.cmdstan/cmdstan-2.35.0")

library(tidyverse)
library(R.utils)
library(wbCorr)
library(readxl)
library(kableExtra)
library(brms)
library(bayesplot)
library(beepr)
library(DHARMa)



source(file.path('Functions', 'ReportModels.R'))
source(file.path('Functions', 'PrettyTables.R'))
source(file.path('Functions', 'ReportMeasures.R'))
source(file.path('Functions', 'PrepareData.R'))
system("shutdown /a")
## [1] 1116
# Set options for analysis
use_mi = FALSE
shutdown = FALSE
report_ordinal = FALSE

options(
  dplyr.print_max = 100, 
  brms.backend = 'cmdstan',
  brms.file_refit = ifelse(use_mi, 'never', 'on_change'),
  brms.file_refit = 'on_change',
  #brms.file_refit = 'always',
  error = function() beepr::beep(sound = 5),
  es.use_symbols = TRUE
)


####################### Model parameters #######################

iterations = 5000 # 10'000 per chain to achieve 40'000
warmup = 2000

# NO AR!!!
#corstr = 'ar'
#corstr = 'cosy_couple'
#corstr = 'cosy_couple:user'


################################################################

suffix = as.character(iterations)
df <- openxlsx::read.xlsx(file.path('long.xlsx'))
df_original <- df

df_double <- prepare_data(df, recode_pushing = TRUE, use_mi = use_mi)[[1]]

Constructing scales Re-coding pusing reshaping data (4field) centering data within and between

Modelling

# For indistinguishable Dyads
model_rows_fixed <- c(
    'Intercept', 
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb'
  )


model_rows_fixed_ordinal <- c(
  model_rows_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rows_fixed[2:length(model_rows_fixed)]
)

model_rows_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)

model_rows_random_ordinal <- c(model_rows_random,'disc')
# For indistinguishable Dyads
model_rownames_fixed <- c(
    'Intercept', 
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'Daily received persuasion target -> target', 
    'Daily received persuasion target -> agent', 
    'Daily received pressure target -> target', 
    'Daily received pressure target -> agent', 
    'Daily received pushing target -> target', 
    'Daily received pushing target -> agent', 
    'Day', 
    'Daily weartime',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'Mean received persuasion target -> target',
    'Mean received persuasion target -> agent',
    'Mean received pressure target -> target',
    'Mean received pressure target -> agent',
    'Mean received pushing target -> target',
    'Mean received pushing target -> agent',
    'Mean weartime'
  )


model_rownames_fixed_ordinal <- c(
  model_rownames_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed[2:length(model_rownames_fixed)]
)

model_rownames_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(Daily received persuasion target -> target)', 
  'sd(Daily received persuasion target -> agent)', 
  'sd(Daily received pressure target -> target)', 
  'sd(Daily received pressure target -> agent)', 
  'sd(Daily received pushing target -> target)', 
  'sd(Daily received pushing target -> agent)', 
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]',
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)

model_rownames_random_ordinal <- c(model_rownames_random,'disc')
# For indistinguishable Dyads
model_rownames_fixed <- c(
    "Intercept", 
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime"
  )


model_rownames_fixed_ordinal <- c(
  model_rownames_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed[2:length(model_rownames_fixed)]
)

model_rownames_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)

model_rownames_random_ordinal <- c(model_rownames_random,'disc')
rows_to_pack <- list(
  "Within-Person Effects" = c(2,9),
  "Between-Person Effects" = c(10,16),
  "Random Effects" = c(17, 23), 
  "Additional Parameters" = c(24,30)
  )


rows_to_pack_ordinal <- list(
  "Intercepts" = c(1,6),
  "Within-Person Effects" = c(2+5,9+5),
  "Between-Person Effects" = c(10+5,16+5),
  "Random Effects" = c(17+5, 23+5), 
  "Additional Parameters" = c(24+5,30+6)
  )

HURDLE MODELS

# For indistinguishable Dyads
model_rows_fixed_hu <- c(
    'Intercept', 
    'hu_Intercept',
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb',
    
    # HURDLE MODEL
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'hu_persuasion_self_cw', 
    'hu_persuasion_partner_cw', 
    'hu_pressure_self_cw', 
    'hu_pressure_partner_cw', 
    'hu_pushing_self_cw', 
    'hu_pushing_partner_cw', 
    'hu_day', 
    'hu_weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'hu_persuasion_self_cb',
    'hu_persuasion_partner_cb',
    'hu_pressure_self_cb',
    'hu_pressure_partner_cb',
    'hu_pushing_self_cb',
    'hu_pushing_partner_cb',
    'hu_weartime_self_cb'
  )

model_rows_random_hu <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(hu_Intercept)',
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # HURDLE
  'sd(hu_persuasion_self_cw)',
  'sd(hu_persuasion_partner_cw)',
  'sd(hu_pressure_self_cw)',
  'sd(hu_pressure_partner_cw)',
  'sd(hu_pushing_self_cw)',
  'sd(hu_pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)
# For indistinguishable Dyads
model_rownames_fixed_hu <- c(
    "Intercept", 
    "Hurdle Intercept",
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime",
    
    # HURDLE
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Hu Daily persuasion experienced", 
    "Hu Daily persuasion utilized (partner's view)", # OR partner received
    "Hu Daily pressure experienced", 
    "Hu Daily pressure utilized (partner's view)", 
    "Hu Daily pushing experienced", 
    "Hu Daily pushing utilized (partner's view)", 
    "Hu Day", 
    "Hu Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Hu Mean persuasion experienced", 
    "Hu Mean persuasion utilized (partner's view)", 
    "Hu Mean pressure experienced", 
    "Hu Mean pressure utilized (partner's view)", 
    "Hu Mean pushing experienced", 
    "Hu Mean pushing utilized (partner's view)", 
    "Hu Mean weartime"
  )


model_rownames_random_hu <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(Hurdle Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  
  # Hurdle
  "sd(Hu Daily persuasion experienced)", 
  "sd(Hu Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Hu Daily pressure experienced)", 
  "sd(Hu Daily pressure utilized (partner's view))", 
  "sd(Hu Daily pushing experienced)", 
  "sd(Hu Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'ar[1]', 
  'ma[1]', 
  'cosy',
  'nu',
  'shape',
  'sderr',
  'sigma'
)
rows_to_pack_hu <- list(
  "Conditional Within-Person Effects" = c(3,10),
  "Conditional Between-Person Effects" = c(11,17),
  
  "Hurdle Within-Person Effects" = c(18,25),
  "Hurdle Between-Person Effects" = c(26,32),
  
  "Random Effects" = c(33, 46), 
  "Additional Parameters" = c(47,53)
  )

Subjective MVPA

range(df_double$pa_sub, na.rm = T) 
## [1]   0 720
hist(df_double$pa_sub, breaks = 40) 

hist(log(df_double$pa_sub+00000000001), breaks = 40)

Hurdle Lognormal Model

formula <- bf(
  pa_sub ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID),
  
  hu = ~ persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 1.5)", class = "b")
  , brms::set_prior("normal(0, 2.5", class = "b", dpar = "hu")
  , brms::set_prior("normal(0, 50)", class = "Intercept") # for non-zero PA
  , brms::set_prior("normal(6, 2.5)", class = "Intercept", dpar = 'hu') # hurdle part
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  
  #, brms::set_prior("normal(10, 10", class = "shape")
  
  #, brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 2)", class = "sderr", lb = 0)
  
  #, autocor_prior
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal_NOAR", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(pa_sub, log_pp_check = TRUE)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 1 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 12000 by 3736 log-likelihood matrix.
## 
##          Estimate    SE
## elpd_loo -10490.4 166.5
## p_loo       181.4   6.2
## looic     20980.8 332.9
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.6, 1.7]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     3735  100.0%  605     
##    (0.7, 1]   (bad)         1    0.0%  <NA>    
##    (1, Inf)   (very bad)    0    0.0%  <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(pa_sub, integer = TRUE, outliers_type = 'bootstrap')

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 9, observations = 3736, p-value = 0.34
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.0005353319 0.0032119914
## sample estimates:
## outlier frequency (expected: 0.00164346895074946 ) 
##                                        0.002408994
#shinystan::launch_shinystan(pa_sub)

In this instance, the warning about max treedepth is a false positive. We have set treedepth to 15, and when we check with shinystan, we see that treedepth is consistently between 10 and 14.

summarize_brms(
  pa_sub, 
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_hu
)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
## Warning in summarize_brms(pa_sub, model_rows_fixed = model_rows_fixed_hu, :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercept 47.90* 42.22 54.30 1.001 3264.34 5142.29 1.000 >100 Overwhelming Evidence
Hurdle Intercept 1.21 0.87 1.69 1.004 2283.61 4134.98 0.868 0.010 Very Strong Evidence for Null
Conditional Within-Person Effects
Daily persuasion experienced 1.03 0.97 1.08 1.000 4843.30 7119.87 0.828 0.028 Very Strong Evidence for Null
Daily persuasion utilized (partner’s view) 1.03 0.98 1.08 1.000 7049.70 8420.12 0.899 0.036 Strong Evidence for Null
Daily pressure experienced 0.89* 0.80 0.99 1.000 10401.96 7267.03 0.984 0.485 Weak Evidence for Null
Daily pressure utilized (partner’s view) 0.94 0.86 1.03 1.000 9300.98 8275.12 0.915 0.075 Strong Evidence for Null
Daily pushing experienced 1.03 0.96 1.10 1.001 7999.07 8660.16 0.775 0.032 Strong Evidence for Null
Daily pushing utilized (partner’s view) 0.99 0.93 1.05 1.001 8640.96 8974.16 0.618 0.022 Very Strong Evidence for Null
Day 1.01 0.89 1.13 1.000 15082.98 9601.07 0.551 0.041 Strong Evidence for Null
Daily weartime NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.01 0.74 1.38 1.000 2732.50 4099.77 0.533 0.102 Moderate Evidence for Null
Mean persuasion utilized (partner’s view) 0.98 0.72 1.33 1.001 2757.84 4589.86 0.543 0.104 Moderate Evidence for Null
Mean pressure experienced 1.15 0.80 1.65 1.001 4069.38 5900.07 0.772 0.162 Moderate Evidence for Null
Mean pressure utilized (partner’s view) 0.89 0.62 1.29 1.001 4198.01 6525.04 0.744 0.150 Moderate Evidence for Null
Mean pushing experienced 1.32 0.84 2.07 1.000 3876.54 6233.24 0.891 0.325 Weak Evidence for Null
Mean pushing utilized (partner’s view) 1.40 0.88 2.21 1.001 3732.45 6367.29 0.925 0.451 Weak Evidence for Null
Mean weartime NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 0.65* 0.57 0.73 1.000 7700.05 7679.91 1.000 >100 Overwhelming Evidence
Hu Daily persuasion utilized (partner’s view) 0.75* 0.67 0.84 1.001 9038.94 9018.25 1.000 >100 Overwhelming Evidence
Hu Daily pressure experienced 1.23 0.88 1.71 1.000 8879.51 8110.20 0.898 0.154 Moderate Evidence for Null
Hu Daily pressure utilized (partner’s view) 0.66* 0.42 0.95 1.000 8465.74 6511.93 0.988 0.899 Weak Evidence for Null
Hu Daily pushing experienced 0.58* 0.40 0.78 1.000 7544.89 7513.38 1.000 24.845 Strong Evidence
Hu Daily pushing utilized (partner’s view) 0.54* 0.41 0.69 1.000 7450.70 6837.16 1.000 >100 Overwhelming Evidence
Hu Day 1.09 0.85 1.42 1.000 15446.28 9424.26 0.750 0.065 Strong Evidence for Null
Hu Daily weartime NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 0.83 0.36 1.87 1.001 1923.05 3521.20 0.680 0.184 Moderate Evidence for Null
Hu Mean persuasion utilized (partner’s view) 0.83 0.36 1.89 1.001 1908.01 3480.06 0.671 0.183 Moderate Evidence for Null
Hu Mean pressure experienced 3.31* 1.36 8.35 1.000 2952.14 5126.23 0.995 6.166 Moderate Evidence
Hu Mean pressure utilized (partner’s view) 1.79 0.74 4.48 1.000 3052.79 5306.39 0.905 0.427 Weak Evidence for Null
Hu Mean pushing experienced 0.35 0.11 1.11 1.002 2801.73 4512.72 0.962 1.144 Weak Evidence
Hu Mean pushing utilized (partner’s view) 0.34 0.11 1.10 1.002 2692.83 4957.40 0.964 1.280 Weak Evidence
Hu Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.32 0.24 0.42 1.00 3342.49 5207.28 NA NA NA
sd(Hurdle Intercept) 0.90 0.69 1.17 1.00 3026.08 5383.96 NA NA NA
sd(Daily persuasion experienced) 0.12 0.08 0.17 1.00 5467.59 7369.26 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.09 0.05 0.13 1.00 6402.38 7446.13 NA NA NA
sd(Daily pressure experienced) 0.08 0.00 0.24 1.00 4554.51 5172.98 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.07 0.00 0.19 1.00 5546.07 5814.39 NA NA NA
sd(Daily pushing experienced) 0.11 0.04 0.19 1.00 5412.48 4100.25 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.09 0.02 0.17 1.00 3656.04 2742.82 NA NA NA
sd(Hu Daily persuasion experienced) 0.18 0.02 0.34 1.00 3562.02 3269.27 NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) 0.17 0.03 0.33 1.00 4124.75 3957.26 NA NA NA
sd(Hu Daily pressure experienced) 0.31 0.01 0.89 1.00 3857.92 3643.66 NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) 0.34 0.01 0.99 1.00 4446.62 5916.04 NA NA NA
sd(Hu Daily pushing experienced) 0.64 0.32 1.07 1.00 5132.64 7726.03 NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) 0.32 0.05 0.64 1.00 4284.77 4067.55 NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma 0.68 0.66 0.71 1.00 15422.04 8395.04 NA NA NA 
mcmc_plot(pa_sub,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

mcmc_plot(pa_sub,
          variable = c(
            'b_hu_persuasion_self_cw',
            'b_hu_persuasion_partner_cw',
            'b_hu_pressure_self_cw',
            'b_hu_pressure_partner_cw',
            'b_hu_pushing_self_cw',
            'b_hu_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95
          )

plot(
  bayestestR::p_direction(pa_sub),
  priors = TRUE
) + 
  #coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in get_priors.brmsfit(model, effects = effects, component = component,
## : NAs introduced by coercion
## Warning in get_priors.brmsfit(model, effects = effects, component = component,
## : NAs introduced by coercion
## Warning: `b_hu_persuasion_self_cw`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_persuasion_partner_cw`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_pressure_self_cw`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_pressure_partner_cw`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_pushing_self_cw`, or one of its groups specified in `by`, is empty
##   and has no density information.
## Warning: `b_hu_pushing_partner_cw`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_persuasion_self_cb`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_persuasion_partner_cb`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_pressure_self_cb`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_pressure_partner_cb`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_pushing_self_cb`, or one of its groups specified in `by`, is empty
##   and has no density information.
## Warning: `b_hu_pushing_partner_cb`, or one of its groups specified in `by`, is
##   empty and has no density information.
## Warning: `b_hu_day`, or one of its groups specified in `by`, is empty and has no
##   density information.
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

conditional_spaghetti(
  pa_sub, 
  effects = c('persuasion_self_cw',
              'persuasion_partner_cw',
              'pressure_self_cw',
              'pressure_partner_cw',
              'pushing_self_cw',
              'pushing_partner_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

\(persuasion_self_cw_count <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-4.png" width="2400" />\)persuasion_self_cw_hurdle \(persuasion_self_cw_combined <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-6.png" width="2400" />\)persuasion_partner_cw_count \(persuasion_partner_cw_hurdle <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-8.png" width="2400" />\)persuasion_partner_cw_combined \(pressure_self_cw_count <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-10.png" width="2400" />\)pressure_self_cw_hurdle \(pressure_self_cw_combined <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-12.png" width="2400" />\)pressure_partner_cw_count \(pressure_partner_cw_hurdle <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-14.png" width="2400" />\)pressure_partner_cw_combined \(pushing_self_cw_count <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-16.png" width="2400" />\)pushing_self_cw_hurdle \(pushing_self_cw_combined <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-18.png" width="2400" />\)pushing_partner_cw_count \(pushing_partner_cw_hurdle <img src="01_FinalModels_files/figure-html/report_pa_sub_lognormal-20.png" width="2400" />\)pushing_partner_cw_combined 

marginaleffects::avg_predictions(pa_sub)
## 
##  Estimate 2.5 % 97.5 %
##      30.9  29.4   32.5
## 
## Type:  response 
## Columns: estimate, conf.low, conf.high

Additionally, as a sensitivity analysis, we estimate the two part of the models separately.

Two separate models

Modelling 0/1

The zero vs. 1 modelling part also has high pareto-k values, but reaches the same conslucsions as the hurdle model. We tried further simplifying by removing the residual AR1 correlation structure, which led to a model with good pareto-k values, still arriving at the same conslusion as the original hurdle model:

df_double$pa_sub_zero <- as.factor(ifelse(df_double$pa_sub > 0, 1, 0))

formula <- bf(
  pa_sub_zero ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
) 

prior1 <- c(
    set_prior("normal(0, 1.5)", class = "b")
  , brms::set_prior("normal(0, 5)", class = "Intercept") # for non-zero PA

  , brms::set_prior("normal(0, 1)", class = "sd", group = "coupleID")

  #, brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 1.5)", class = "sderr", lb = 0)
  
  #, autocor_prior
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )


pa_sub_zero_model <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  #family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_zero_part_NOAR", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(pa_sub_zero_model, log_pp_check = FALSE)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 2 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 12000 by 3736 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo  -2149.7 27.8
## p_loo        93.4  4.1
## looic      4299.5 55.6
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.5, 2.3]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     3734  99.9%   355     
##    (0.7, 1]   (bad)         2   0.1%   <NA>    
##    (1, Inf)   (very bad)    0   0.0%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(pa_sub_zero_model, integer = TRUE, outliers_type = 'bootstrap')

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 0, observations = 3736, p-value = 1
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0 0
## sample estimates:
## outlier frequency (expected: 0 ) 
##                                0
summarize_brms(
  pa_sub_zero_model, 
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_hu
)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
OR l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercept 0.84 0.61 1.15 1.001 1841.97 3973.90 0.863 0.059 Strong Evidence for Null
Hurdle Intercept NA NA NA NA NA NA NA NA NA
Conditional Within-Person Effects
Daily persuasion experienced 1.53* 1.36 1.75 1.001 8062.00 8582.55 1.000 >100 Overwhelming Evidence
Daily persuasion utilized (partner’s view) 1.33* 1.19 1.50 1.000 9007.75 8833.49 1.000 >100 Overwhelming Evidence
Daily pressure experienced 0.82 0.59 1.11 1.000 10666.20 8323.82 0.906 0.242 Moderate Evidence for Null
Daily pressure utilized (partner’s view) 1.48* 1.05 2.25 1.000 10016.51 6740.30 0.987 1.589 Weak Evidence
Daily pushing experienced 1.72* 1.28 2.43 1.000 7256.09 7373.66 1.000 37.284 Very Strong Evidence
Daily pushing utilized (partner’s view) 1.83* 1.46 2.40 1.001 8929.39 8028.12 1.000 >100 Overwhelming Evidence
Day 0.92 0.71 1.19 1.000 18939.10 8954.86 0.739 0.102 Moderate Evidence for Null
Daily weartime NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.17 0.56 2.48 1.002 1891.71 4153.21 0.667 0.270 Moderate Evidence for Null
Mean persuasion utilized (partner’s view) 1.17 0.57 2.44 1.002 1970.51 4133.98 0.666 0.269 Moderate Evidence for Null
Mean pressure experienced 0.33* 0.15 0.75 1.001 2989.43 5922.42 0.996 9.942 Moderate Evidence
Mean pressure utilized (partner’s view) 0.59 0.26 1.34 1.002 3409.72 5951.96 0.900 0.616 Weak Evidence for Null
Mean pushing experienced 2.45 0.84 7.10 1.001 3185.19 5355.11 0.953 1.480 Weak Evidence
Mean pushing utilized (partner’s view) 2.56 0.88 7.30 1.001 3160.95 5396.08 0.958 1.614 Weak Evidence
Mean weartime NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Hu Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Daily pressure experienced NA NA NA NA NA NA NA NA NA
Hu Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Daily pushing experienced NA NA NA NA NA NA NA NA NA
Hu Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Day NA NA NA NA NA NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Hu Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean pressure experienced NA NA NA NA NA NA NA NA NA
Hu Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean pushing experienced NA NA NA NA NA NA NA NA NA
Hu Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.89 0.68 1.14 1.00 3199.26 5243.78 NA NA NA
sd(Hurdle Intercept) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion experienced) 0.18 0.02 0.34 1.00 3164.76 2298.66 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.17 0.02 0.33 1.00 3506.22 3289.23 NA NA NA
sd(Daily pressure experienced) 0.28 0.01 0.78 1.00 4188.59 5081.55 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.31 0.01 0.90 1.00 3962.03 5276.36 NA NA NA
sd(Daily pushing experienced) 0.62 0.32 1.01 1.00 5746.41 6689.15 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.31 0.05 0.63 1.00 4373.04 3371.61 NA NA NA
sd(Hu Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma NA NA NA NA NA NA NA NA NA 
brms::mcmc_plot(pa_sub_zero_model,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

plot(
  bayestestR::p_direction(pa_sub_zero_model),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(pa_sub_zero_model)
) + theme_bw()
## Possible multicollinearity between b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.81). This might lead to inappropriate
##   results. See 'Details' in '?rope'.

conditional_spaghetti(
  pa_sub_zero_model, 
  effects = c('persuasion_self_cw',
              'persuasion_partner_cw',
              'pressure_partner_cw',
              'pushing_self_cw',
              'pushing_partner_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

\(persuasion_self_cw <img src="01_FinalModels_files/figure-html/report_pa_sub_separate_inactive-4.png" width="2400" />\)persuasion_partner_cw \(pressure_partner_cw <img src="01_FinalModels_files/figure-html/report_pa_sub_separate_inactive-6.png" width="2400" />\)pushing_self_cw $pushing_partner_cw 

marginaleffects::avg_predictions(pa_sub_zero_model)
## 
##  Estimate 2.5 % 97.5 %
##     0.447 0.434  0.461
## 
## Type:  response 
## Columns: estimate, conf.low, conf.high

Modelling active Days

# Only include active days
df_double$pa_sub_non_zero <- ifelse(df_double$pa_sub > 0, df_double$pa_sub, NA)

formula <- bf(
  log(pa_sub_non_zero) ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  ##, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 50)", class = "Intercept") # for non-zero PA

  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  
  , brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 2)", class = "sderr", lb = 0)
  
  #, autocor_prior
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )


pa_sub_active_model <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_active_part_NOAR", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(pa_sub_active_model, log_pp_check = TRUE)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
## Computed from 12000 by 1672 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo  -1785.5 35.3
## p_loo        85.9  4.1
## looic      3571.1 70.6
## ------
## MCSE of elpd_loo is 0.1.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.5, 2.0]).
## 
## All Pareto k estimates are good (k < 0.7).
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(pa_sub_active_model, integer = TRUE, outliers_type = 'bootstrap')
## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 8, observations = 1672, p-value = 0.26
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.0005980861 0.0068929426
## sample estimates:
## outlier frequency (expected: 0.00319377990430622 ) 
##                                        0.004784689
summarize_brms(
  pa_sub_active_model, 
  model_rows_fixed = model_rows_fixed_hu,
  model_rows_random = model_rows_random_hu,
  model_rownames_fixed = model_rownames_fixed_hu,
  model_rownames_random = model_rownames_random_hu,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_hu
)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
## Warning in summarize_brms(pa_sub_active_model, model_rows_fixed =
## model_rows_fixed_hu, : Coefficients were exponentiated. Double check if this
## was intended.
exp(Est.) l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercept 47.85* 42.22 54.17 1.001 4041.60 6381.47 1.000 >100 Overwhelming Evidence
Hurdle Intercept NA NA NA NA NA NA NA NA NA
Conditional Within-Person Effects
Daily persuasion experienced 1.03 0.97 1.08 1.000 6597.52 8346.37 0.833 0.018 Very Strong Evidence for Null
Daily persuasion utilized (partner’s view) 1.03 0.99 1.08 1.000 8955.61 9602.22 0.904 0.019 Very Strong Evidence for Null
Daily pressure experienced 0.89* 0.80 0.99 1.000 14496.32 8806.27 0.987 0.287 Moderate Evidence for Null
Daily pressure utilized (partner’s view) 0.94 0.85 1.03 1.000 13817.24 8187.83 0.908 0.044 Strong Evidence for Null
Daily pushing experienced 1.03 0.96 1.10 1.000 11709.14 9409.07 0.775 0.018 Very Strong Evidence for Null
Daily pushing utilized (partner’s view) 0.99 0.93 1.05 1.000 11995.74 10105.19 0.633 0.014 Very Strong Evidence for Null
Day 1.01 0.89 1.14 1.000 20035.93 9281.40 0.545 0.025 Very Strong Evidence for Null
Daily weartime NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.00 0.73 1.37 1.001 3231.81 4989.35 0.510 0.063 Strong Evidence for Null
Mean persuasion utilized (partner’s view) 0.98 0.72 1.33 1.001 3286.72 5489.83 0.560 0.061 Strong Evidence for Null
Mean pressure experienced 1.15 0.80 1.64 1.000 5002.27 7220.72 0.778 0.097 Strong Evidence for Null
Mean pressure utilized (partner’s view) 0.89 0.61 1.28 1.000 5018.40 7058.04 0.732 0.088 Strong Evidence for Null
Mean pushing experienced 1.34 0.84 2.11 1.000 4624.61 7355.15 0.896 0.206 Moderate Evidence for Null
Mean pushing utilized (partner’s view) 1.42 0.89 2.26 1.000 4706.03 7431.81 0.929 0.289 Moderate Evidence for Null
Mean weartime NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced NA NA NA NA NA NA NA NA NA
Hu Daily persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Daily pressure experienced NA NA NA NA NA NA NA NA NA
Hu Daily pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Daily pushing experienced NA NA NA NA NA NA NA NA NA
Hu Daily pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Day NA NA NA NA NA NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced NA NA NA NA NA NA NA NA NA
Hu Mean persuasion utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean pressure experienced NA NA NA NA NA NA NA NA NA
Hu Mean pressure utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean pushing experienced NA NA NA NA NA NA NA NA NA
Hu Mean pushing utilized (partner’s view) NA NA NA NA NA NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.32 0.25 0.42 1.00 3621.92 6455.57 NA NA NA
sd(Hurdle Intercept) NA NA NA NA NA NA NA NA NA
sd(Daily persuasion experienced) 0.12 0.08 0.17 1.00 7194.08 9303.03 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.09 0.05 0.13 1.00 8556.16 8646.99 NA NA NA
sd(Daily pressure experienced) 0.08 0.00 0.23 1.00 5748.83 6157.17 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.07 0.00 0.19 1.00 6534.51 6674.27 NA NA NA
sd(Daily pushing experienced) 0.11 0.04 0.19 1.00 5752.73 4564.11 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.09 0.02 0.17 1.00 5084.36 3390.77 NA NA NA
sd(Hu Daily persuasion experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing experienced) NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) NA NA NA NA NA NA NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma 0.68 0.66 0.71 1.00 19220.03 8549.91 NA NA NA 
mcmc_plot(pa_sub_active_model,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

plot(
  bayestestR::p_direction(pa_sub_active_model),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(pa_sub_active_model)
) + theme_bw()
## Possible multicollinearity between b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.77). This might lead to inappropriate
##   results. See 'Details' in '?rope'.

conditional_spaghetti(
  pa_sub_active_model, 
  effects = c('pressure_self_cw'),
  transform_fn = function(x) {exp(x)},
  group_var = 'coupleID',
  plot_full_range = TRUE
)

$pressure_self_cw 

marginaleffects::avg_predictions(pa_sub_active_model)
## 
##  Estimate 2.5 % 97.5 %
##      3.92  3.89   3.95
## 
## Type:  response 
## Columns: estimate, conf.low, conf.high

Device Based MVPA

range(df_double$pa_obj, na.rm = T) 
## [1]   5.75 971.25
hist(df_double$pa_obj, breaks = 50)

df_double$pa_obj_log <- log(df_double$pa_obj)

hist(df_double$pa_obj_log, breaks = 50)

log gaussian

formula <- bf(
  log(pa_obj) ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day + weartime_self_cw + weartime_self_cb +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 50)", class = "Intercept") # for non-zero PA

  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  
  , brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 2)", class = "sderr", lb = 0)
  
  #, autocor_prior
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian_NOAR", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(pa_obj_log, log_pp_check = TRUE)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
## Computed from 12000 by 3337 log-likelihood matrix.
## 
##          Estimate    SE
## elpd_loo  -2936.5  56.3
## p_loo        92.1   4.5
## looic      5872.9 112.6
## ------
## MCSE of elpd_loo is 0.1.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.4, 2.3]).
## 
## All Pareto k estimates are good (k < 0.7).
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(pa_obj_log, integer = TRUE, outliers_type = 'bootstrap')
## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 26, observations = 3337, p-value < 2.2e-16
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.001341025 0.005251723
## sample estimates:
## outlier frequency (expected: 0.00302966736589751 ) 
##                                        0.007791429
summarize_brms(
  pa_obj_log, 
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
## Warning in summarize_brms(pa_obj_log, model_rows_fixed = model_rows_fixed, :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercept 117.41* 105.48 130.41 1.002 1936.20 4337.49 1.000 >100 Overwhelming Evidence
Within-Person Effects
Daily persuasion experienced 1.03 1.00 1.06 1.000 8443.32 9440.38 0.966 0.034 Strong Evidence for Null
Daily persuasion utilized (partner’s view) 1.02 0.99 1.05 1.001 10110.26 9523.00 0.888 0.014 Very Strong Evidence for Null
Daily pressure experienced 0.94 0.88 1.01 1.000 13123.29 8882.26 0.960 0.073 Strong Evidence for Null
Daily pressure utilized (partner’s view) 0.98 0.92 1.05 1.000 14150.05 9457.65 0.714 0.015 Very Strong Evidence for Null
Daily pushing experienced 1.03 0.98 1.08 1.000 11136.01 8673.71 0.900 0.024 Very Strong Evidence for Null
Daily pushing utilized (partner’s view) 1.02 0.97 1.06 1.000 15288.78 9998.46 0.771 0.011 Very Strong Evidence for Null
Day 0.97 0.91 1.04 1.000 20771.69 8890.89 0.785 0.019 Very Strong Evidence for Null
Daily weartime 1.00* 1.00 1.00 1.000 12258.92 8147.17 1.000 0.067 Strong Evidence for Null
Between-Person Effects
Mean persuasion experienced 1.10 0.82 1.46 1.003 1635.21 3066.39 0.748 0.074 Strong Evidence for Null
Mean persuasion utilized (partner’s view) 0.98 0.73 1.30 1.003 1664.34 3415.83 0.559 0.059 Strong Evidence for Null
Mean pressure experienced 0.98 0.73 1.31 1.003 2233.86 4590.15 0.560 0.060 Strong Evidence for Null
Mean pressure utilized (partner’s view) 0.97 0.73 1.28 1.002 2141.33 4170.53 0.594 0.060 Strong Evidence for Null
Mean pushing experienced 0.97 0.65 1.45 1.002 2615.53 4583.92 0.553 0.083 Strong Evidence for Null
Mean pushing utilized (partner’s view) 1.25 0.83 1.85 1.002 2546.44 4664.09 0.859 0.146 Moderate Evidence for Null
Mean weartime 1.00 1.00 1.00 1.000 16930.39 10509.93 0.909 0.000 Very Strong Evidence for Null
Random Effects
sd(Intercept) 0.31 0.24 0.40 1.00 2373.00 4540.45 NA NA NA
sd(Daily persuasion experienced) 0.05 0.02 0.08 1.00 6108.14 4438.88 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.06 0.03 0.09 1.00 6116.29 5901.98 NA NA NA
sd(Daily pressure experienced) 0.05 0.00 0.14 1.00 5104.76 5590.89 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.04 0.00 0.12 1.00 6898.74 5731.73 NA NA NA
sd(Daily pushing experienced) 0.07 0.01 0.15 1.00 2745.29 3354.31 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.04 0.00 0.10 1.00 4176.37 5266.68 NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma 0.57 0.56 0.59 1.00 19512.21 8255.47 NA NA NA 
mcmc_plot(pa_obj_log,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

plot(
  bayestestR::p_direction(pa_obj_log),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(pa_obj_log)
) + theme_bw()
## Possible multicollinearity between b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.89), b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.74), b_pressure_partner_cb and
##   b_persuasion_self_cb (r = 0.74), b_pressure_self_cb and
##   b_persuasion_partner_cb (r = 0.72), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.78), b_pushing_partner_cb and
##   b_pushing_self_cb (r = 0.82). This might lead to inappropriate results.
##   See 'Details' in '?rope'.

# Nothing significant, no plots

Affect

range(df_double$aff, na.rm = T) 
## [1] 0 5
hist(df_double$aff, breaks = 15)

Gaussian

formula <- bf(
  aff ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  ,brms::set_prior("normal(0, 20)", class = "Intercept", lb=1, ub=6) # range of the outcome scale
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  
  , brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 2)", class = "sderr", lb = 0)
  
  #, autocor_prior
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss_NOAR", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(mood_gauss, log_pp_check = FALSE)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
## Computed from 12000 by 3736 log-likelihood matrix.
## 
##          Estimate    SE
## elpd_loo  -5185.9  59.2
## p_loo        75.2   3.3
## looic     10371.7 118.5
## ------
## MCSE of elpd_loo is 0.1.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.5, 2.1]).
## 
## All Pareto k estimates are good (k < 0.7).
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(mood_gauss, integer = FALSE)

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 31, observations = 3736, p-value =
## 9.752e-11
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.00564461 0.01175733
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.008297645
summarize_brms(
  mood_gauss, 
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = F) %>%
  print_df(rows_to_pack = rows_to_pack)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
b l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercept 3.70* 3.48 3.91 1.007 1299.69 2626.71 1.000 >100 Overwhelming Evidence
Within-Person Effects
Daily persuasion experienced 0.00 -0.04 0.05 1.000 10096.58 8360.31 0.553 0.004 Very Strong Evidence for Null
Daily persuasion utilized (partner’s view) 0.02 -0.02 0.07 1.001 8776.31 8901.48 0.829 0.007 Very Strong Evidence for Null
Daily pressure experienced -0.04 -0.14 0.07 1.000 11418.22 8197.57 0.767 0.013 Very Strong Evidence for Null
Daily pressure utilized (partner’s view) -0.02 -0.14 0.08 1.000 10703.53 8437.61 0.661 0.012 Very Strong Evidence for Null
Daily pushing experienced 0.02 -0.04 0.09 1.000 11578.47 8985.61 0.768 0.008 Very Strong Evidence for Null
Daily pushing utilized (partner’s view) 0.08* 0.01 0.14 1.000 10076.04 8287.53 0.984 0.078 Strong Evidence for Null
Day 0.26* 0.15 0.37 1.001 17305.49 8858.10 1.000 50.970 Very Strong Evidence
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 0.33 -0.21 0.90 1.003 1324.40 2167.41 0.889 0.117 Moderate Evidence for Null
Mean persuasion utilized (partner’s view) 0.22 -0.32 0.79 1.003 1322.88 2252.89 0.790 0.079 Strong Evidence for Null
Mean pressure experienced -0.31 -0.86 0.23 1.003 1481.77 2769.72 0.874 0.102 Moderate Evidence for Null
Mean pressure utilized (partner’s view) -0.30 -0.84 0.24 1.003 1476.58 2758.08 0.871 0.098 Strong Evidence for Null
Mean pushing experienced 0.22 -0.54 0.99 1.002 1997.17 3840.31 0.714 0.096 Strong Evidence for Null
Mean pushing utilized (partner’s view) 0.35 -0.39 1.12 1.002 2005.43 3822.38 0.825 0.119 Moderate Evidence for Null
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.60 0.47 0.78 1.00 2148.52 4074.04 NA NA NA
sd(Daily persuasion experienced) 0.04 0.00 0.10 1.00 3348.96 5109.75 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.08 0.01 0.13 1.00 3084.49 2127.56 NA NA NA
sd(Daily pressure experienced) 0.08 0.00 0.24 1.00 4964.63 5526.49 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.09 0.00 0.26 1.00 4911.35 5368.11 NA NA NA
sd(Daily pushing experienced) 0.06 0.00 0.14 1.00 3656.74 3408.48 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.07 0.00 0.17 1.00 4184.55 4073.31 NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma 0.96 0.94 0.98 1.00 18080.51 9327.96 NA NA NA 
mcmc_plot(mood_gauss,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

plot(
  bayestestR::p_direction(mood_gauss),
  priors = TRUE
)  + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(mood_gauss)
) + theme_bw()
## Possible multicollinearity between b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.82), b_pressure_partner_cb and
##   b_persuasion_self_cb (r = 0.8), b_pressure_self_cb and
##   b_persuasion_partner_cb (r = 0.8), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.82), b_pressure_partner_cb and
##   b_pressure_self_cb (r = 0.76), b_pushing_partner_cb and
##   b_pushing_self_cb (r = 0.89). This might lead to inappropriate results.
##   See 'Details' in '?rope'.

conditional_spaghetti(
  mood_gauss, 
  effects = c('pushing_partner_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

$pushing_partner_cw 

marginaleffects::avg_predictions(mood_gauss)
## 
##  Estimate 2.5 % 97.5 %
##      3.82  3.79   3.85
## 
## Type:  response 
## Columns: estimate, conf.low, conf.high

Reactance

range(df_double$reactance, na.rm = T) 
## [1] 0 5
hist(df_double$reactance, breaks = 7) 

hist(log(df_double$reactance+0.1), breaks = 10)

Ordinal

df_double$reactance_ordinal <- factor(df_double$reactance,
                                      levels = 0:5, 
                                      ordered = TRUE)

formula <- bf(
  reactance_ordinal ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  #, brms::set_prior("normal(0, 10)", class = "Intercept", lb=0, ub=5) # range of the outcome scale
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  
  #, brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 2)", class = "sderr", lb = 0)
  
  #, autocor_prior
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = cumulative() # HURDLE_CUMULATIVE
#  )


#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

reactance_ordinal <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::cumulative(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777
  , file = file.path("models_cache_brms", paste0("reactance_ordinal_NOARNOAR_", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(reactance_ordinal)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 6 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 12000 by 756 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo   -681.9 31.9
## p_loo        73.4  5.4
## looic      1363.8 63.8
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.4, 1.7]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     750   99.2%   445     
##    (0.7, 1]   (bad)        6    0.8%   <NA>    
##    (1, Inf)   (very bad)   0    0.0%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(reactance_ordinal, outliers_type = 'bootstrap')

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 2, observations = 756, p-value = 0.02
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.000000000 0.001322751
## sample estimates:
## outlier frequency (expected: 0.000291005291005291 ) 
##                                         0.002645503
summarize_brms(
  reactance_ordinal, 
  model_rows_fixed = model_rows_fixed_ordinal,
  model_rows_random = model_rows_random_ordinal,
  model_rownames_fixed = model_rownames_fixed_ordinal,
  model_rownames_random = model_rownames_random_ordinal,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_ordinal)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
OR l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercepts
Intercept NA NA NA NA NA NA NA NA NA
Intercept[1] 3.85* 2.33 6.45 1.000 8597.51 9211.52 1.000 >100 Overwhelming Evidence
Intercept[2] 8.35* 4.95 14.45 1.000 8834.13 8761.87 1.000 >100 Overwhelming Evidence
Intercept[3] 23.24* 13.13 42.31 1.000 9373.11 9135.25 1.000 >100 Overwhelming Evidence
Intercept[4] 101.58* 52.10 209.10 1.000 10595.59 9691.96 1.000 >100 Overwhelming Evidence
Intercept[5] 3488.39* 1077.21 13336.21 1.001 13289.90 8614.35 1.000 >100 Overwhelming Evidence
Within-Person Effects
Daily persuasion experienced 0.85* 0.71 0.99 1.000 9896.79 7775.98 0.980 0.249 Moderate Evidence for Null
Daily persuasion utilized (partner’s view) 1.03 0.84 1.24 1.000 9229.64 7954.31 0.607 0.041 Strong Evidence for Null
Daily pressure experienced 1.84* 1.18 2.68 1.001 5608.70 6574.26 0.994 3.418 Moderate Evidence
Daily pressure utilized (partner’s view) 1.22 0.71 2.06 1.001 7420.47 6712.38 0.808 0.153 Moderate Evidence for Null
Daily pushing experienced 1.17 0.97 1.43 1.001 7757.54 7485.59 0.946 0.151 Moderate Evidence for Null
Daily pushing utilized (partner’s view) 0.91 0.71 1.17 1.000 9792.98 8248.45 0.770 0.064 Strong Evidence for Null
Day 1.47 0.75 2.89 1.000 13856.59 9556.08 0.871 0.255 Moderate Evidence for Null
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.12 0.40 3.11 1.000 4004.78 6205.61 0.586 0.209 Moderate Evidence for Null
Mean persuasion utilized (partner’s view) 1.38 0.45 4.31 1.000 4278.63 6604.64 0.711 0.252 Moderate Evidence for Null
Mean pressure experienced 3.51* 1.18 10.71 1.000 4792.44 6629.78 0.990 3.171 Moderate Evidence
Mean pressure utilized (partner’s view) 1.17 0.37 3.66 1.000 4762.77 6973.89 0.612 0.233 Moderate Evidence for Null
Mean pushing experienced 1.23 0.28 5.62 1.000 5266.91 7722.24 0.605 0.312 Weak Evidence for Null
Mean pushing utilized (partner’s view) 0.11* 0.02 0.64 1.000 7169.66 7894.75 0.992 7.489 Moderate Evidence
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.82 0.47 1.25 1.00 4256.66 7026.55 NA NA NA
sd(Daily persuasion experienced) 0.18 0.01 0.43 1.00 1818.16 4469.30 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.22 0.01 0.51 1.00 3273.96 4646.21 NA NA NA
sd(Daily pressure experienced) 0.56 0.09 1.14 1.00 2689.26 2386.75 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.50 0.02 1.55 1.00 2929.20 4951.16 NA NA NA
sd(Daily pushing experienced) 0.22 0.01 0.50 1.00 3370.38 4130.06 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.18 0.01 0.57 1.00 4830.97 5284.06 NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma NA NA NA NA NA NA NA NA NA
disc 1.00 1.00 1.00 NA NA NA NA NA NA 
mcmc_plot(reactance_ordinal,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

plot(
  bayestestR::p_direction(reactance_ordinal),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(reactance_ordinal)
) + theme_bw()
## Possible multicollinearity between b_Intercept[4] and b_Intercept[2] (r
##   = 0.76), b_Intercept[4] and b_Intercept[3] (r = 0.83),
##   b_pressure_self_cb and b_persuasion_self_cb (r = 0.72),
##   b_pressure_partner_cb and b_persuasion_partner_cb (r = 0.79). This might
##   lead to inappropriate results. See 'Details' in '?rope'.

conditional_spaghetti(
  reactance_ordinal, 
  effects = c('persuasion_self_cw', 'pressure_self_cw')
  , group_var = 'coupleID'
  #, n_groups = 15
  , plot_full_range = T
)

\(persuasion_self_cw <img src="01_FinalModels_files/figure-html/report_reactance_ordinal-4.png" width="2400" />\)pressure_self_cw 

marginaleffects::avg_predictions(reactance_ordinal)
## 
##  Group Estimate   2.5 % 97.5 %
##      0  0.68180 0.65445 0.7078
##      1  0.09430 0.07497 0.1157
##      2  0.08365 0.06601 0.1047
##      3  0.06983 0.05482 0.0875
##      4  0.06396 0.05134 0.0776
##      5  0.00531 0.00172 0.0116
## 
## Type:  response 
## Columns: group, estimate, conf.low, conf.high

Binary

introduce_binary_reactance <- function(data) {
  data$is_reactance <- factor(data$reactance > 0, levels = c(FALSE, TRUE), labels = c(0, 1))
  return(data)
}



df_double <- introduce_binary_reactance(df_double)
if (use_mi) {
  for (i in seq_along(implist)) {
    implist[[i]] <- introduce_binary_reactance(implist[[i]])
  }
}


formula <- bf(
  is_reactance ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 10)", class = "Intercept", lb=0, ub=5) # range of the outcome scale
  
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  
  #, brms::set_prior("cauchy(0, 2)", class = "sigma", lb = 0)
  #, brms::set_prior("cauchy(0, 2)", class = "sderr", lb = 0)
  
  #, autocor_prior
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance_NOARNOAR_", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
check_brms(is_reactance)
## 
## Divergences:
## 0 of 12000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 12000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 55 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.
## 
## Computed from 12000 by 756 log-likelihood matrix.
## 
##          Estimate   SE
## elpd_loo   -363.1 16.0
## p_loo        80.1  6.0
## looic       726.3 32.0
## ------
## MCSE of elpd_loo is NA.
## MCSE and ESS estimates assume MCMC draws (r_eff in [0.4, 1.5]).
## 
## Pareto k diagnostic values:
##                          Count Pct.    Min. ESS
## (-Inf, 0.7]   (good)     701   92.7%   429     
##    (0.7, 1]   (bad)       49    6.5%   <NA>    
##    (1, Inf)   (very bad)   6    0.8%   <NA>    
## See help('pareto-k-diagnostic') for details.
DHARMa.check_brms.all(is_reactance, integer = FALSE)

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 1, observations = 756, p-value = 1
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.0000334886 0.0073476538
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.001322751
summarize_brms(
  is_reactance, 
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
OR l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS p_dir BF BF_Evidence
Intercept 0.29* 0.16 0.50 1.000 13412.80 10197.69 1.000 >100 Overwhelming Evidence
Within-Person Effects
Daily persuasion experienced 0.84 0.69 1.01 1.001 12653.33 8961.78 0.966 0.197 Moderate Evidence for Null
Daily persuasion utilized (partner’s view) 1.13 0.85 1.54 1.000 9878.26 8807.43 0.789 0.076 Strong Evidence for Null
Daily pressure experienced 2.04* 1.03 4.57 1.000 8289.60 6530.31 0.979 1.397 Weak Evidence
Daily pressure utilized (partner’s view) 1.44 0.58 4.04 1.000 8810.48 6987.21 0.799 0.241 Moderate Evidence for Null
Daily pushing experienced 1.28* 1.01 1.64 1.001 12865.46 8547.86 0.979 0.398 Weak Evidence for Null
Daily pushing utilized (partner’s view) 0.89 0.60 1.31 1.000 14162.09 8870.67 0.735 0.095 Strong Evidence for Null
Day 1.64 0.77 3.51 1.000 19648.18 9310.83 0.900 0.359 Weak Evidence for Null
Daily weartime NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.99 0.60 6.80 1.000 6422.05 7842.49 0.870 0.443 Weak Evidence for Null
Mean persuasion utilized (partner’s view) 1.90 0.52 7.15 1.000 6799.76 8735.29 0.829 0.409 Weak Evidence for Null
Mean pressure experienced 17.91* 2.33 159.59 1.000 8975.03 8541.37 0.997 20.141 Strong Evidence
Mean pressure utilized (partner’s view) 2.27 0.24 19.14 1.000 8119.78 8948.63 0.769 0.607 Weak Evidence for Null
Mean pushing experienced 0.83 0.12 6.30 1.000 8730.85 8605.26 0.580 0.408 Weak Evidence for Null
Mean pushing utilized (partner’s view) 0.08* 0.01 0.66 1.000 10022.27 10004.57 0.990 7.272 Moderate Evidence
Mean weartime NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 1.18 0.75 1.74 1.00 5382.54 7890.40 NA NA NA
sd(Daily persuasion experienced) 0.21 0.01 0.51 1.00 2939.60 5289.74 NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.50 0.12 0.97 1.00 4502.23 4854.72 NA NA NA
sd(Daily pressure experienced) 1.12 0.16 2.43 1.00 2982.38 3385.24 NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.97 0.04 2.73 1.00 4249.53 5650.95 NA NA NA
sd(Daily pushing experienced) 0.25 0.01 0.60 1.00 4754.94 5282.30 NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.31 0.01 0.95 1.00 4810.38 6592.98 NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA
sigma NA NA NA NA NA NA NA NA NA 
mcmc_plot(is_reactance,
          variable = c(
            'b_persuasion_self_cw',
            'b_persuasion_partner_cw',
            'b_pressure_self_cw',
            'b_pressure_partner_cw',
            'b_pushing_self_cw',
            'b_pushing_partner_cw'
          ),
          #regex = TRUE,
          type = 'areas',
          prob = 0.95)

plot(
  bayestestR::p_direction(is_reactance),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(is_reactance)
) + theme_bw()

conditional_spaghetti(
  is_reactance, 
  effects = c('pressure_self_cw', 'pushing_self_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

\(pressure_self_cw <img src="01_FinalModels_files/figure-html/report_is_reactance-4.png" width="2400" />\)pushing_self_cw 

marginaleffects::avg_predictions(is_reactance)
## 
##  Estimate 2.5 % 97.5 %
##     0.328 0.303  0.355
## 
## Type:  response 
## Columns: estimate, conf.low, conf.high
hypothesis(is_reactance, "pressure_self_cw > pushing_self_cw")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (pressure_self_cw... > 0     0.47      0.39    -0.15     1.13       8.86
##   Post.Prob Star
## 1       0.9     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.

Report All Models

if (report_ordinal) {
  model_rows_random_final <- model_rows_random_ordinal
  model_rows_fixed_final <- model_rows_fixed_ordinal
  model_rownames_fixed_final <- model_rownames_fixed_ordinal
  model_rownames_random_final <- model_rownames_random_ordinal
  rows_to_pack_final <- rows_to_pack_ordinal
} else {
  model_rows_random_final <- model_rows_random_hu
  model_rows_fixed_final <- model_rows_fixed_hu
  model_rownames_fixed_final <- model_rownames_fixed_hu
  model_rownames_random_final <- model_rownames_random_hu
  rows_to_pack_final <- rows_to_pack_hu
}


bayesfactor = TRUE

all_models <- report_side_by_side(
  pa_sub,
  pa_obj_log,
  mood_gauss,
  reactance_ordinal,
  is_reactance,
  
  bayesfactor = bayesfactor,
  model_rows_random = model_rows_random_final,
  model_rows_fixed = model_rows_fixed_final,
  model_rownames_random = model_rownames_random_final,
  model_rownames_fixed = model_rownames_fixed_final
)

[1] “pa_sub”

## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
## Warning in summarize_brms(model, short_version = TRUE, bayesfactor =
## bayesfactor, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “pa_obj_log”

## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
## Warning: Coefficients were exponentiated. Double check if this was intended.

[1] “mood_gauss”

## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.

[1] “reactance_ordinal”

## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.

[1] “is_reactance”

## Sampling priors, please wait...
## Warning: Bayes factors might not be precise.
##   For precise Bayes factors, sampling at least 40,000 posterior samples is
##   recommended.
# pretty printing
if (bayesfactor) {
  summary_all_models <- all_models %>%
    print_df(rows_to_pack = rows_to_pack_final) %>%
    add_header_above(
      c(" ", "Subjective MVPA Hurdle Lognormal" = 5, 
        "Device-Based MVPA Log (Gaussian)" = 5, 
        "Mood Gaussian" = 5,
        "Reactance Ordinal" = 5,
        "Reactance Dichotome" = 5)
    )

  export_xlsx(
    summary_all_models, 
    rows_to_pack = rows_to_pack_final,
    file.path("Output", "AllModels_experimental_noAR.xlsx"), 
    merge_option = 'both', 
    simplify_2nd_row = TRUE,
    colwidths = c(38, 
                  7.2,13.3,7.2,7.2,24,
                  7.2,13.3,7.2,7.2,24,
                  7.2,13.3,7.2,7.2,24,
                  7.2,13.3,7.2,7.2,24,
                  7.2,13.3,7.2,7.2,24),
    line_above_rows = c(1,2),
    line_below_rows = c(-1)
  )
} else {
  summary_all_models <- all_models %>%
    print_df(rows_to_pack = rows_to_pack_final) %>%
    add_header_above(
      c(" ", "Subjective MVPA Hurdle Lognormal" = 3, 
        "Device-Based MVPA Log (Gaussian)" = 3, 
        "Mood Gaussian" = 3,
        "Reactance Ordinal" = 3,
        "Reactance Dichotome" = 3)
    )

  export_xlsx(
    summary_all_models, 
    rows_to_pack = rows_to_pack_final,
    file.path("Output", "AllModels_experimental_noAR.xlsx"), 
    merge_option = 'both', 
    simplify_2nd_row = TRUE,
    colwidths = c(38, 
                  7.2,13.3,7.2,
                  7.2,13.3,7.2,
                  7.2,13.3,7.2,
                  7.2,13.3,7.2,
                  7.2,13.3,7.2),
    line_above_rows = c(1,2),
    line_below_rows = c(-1)
  )
}
## 
## Attaching package: 'rvest'
## The following object is masked from 'package:readr':
## 
##     guess_encoding
print(summary_all_models)
Subjective MVPA Hurdle Lognormal
Device-Based MVPA Log (Gaussian)
Mood Gaussian
Reactance Ordinal
Reactance Dichotome
exp(Est.) pa_sub 95% CI pa_sub p_dir pa_sub BF pa_sub BF_Evidence pa_sub exp(Est.) pa_obj_log 95% CI pa_obj_log p_dir pa_obj_log BF pa_obj_log BF_Evidence pa_obj_log b mood_gauss 95% CI mood_gauss p_dir mood_gauss BF mood_gauss BF_Evidence mood_gauss OR reactance_ordinal 95% CI reactance_ordinal p_dir reactance_ordinal BF reactance_ordinal BF_Evidence reactance_ordinal OR is_reactance 95% CI is_reactance p_dir is_reactance BF is_reactance BF_Evidence is_reactance
Intercept 47.90* [42.22, 54.30] 1.000 >100 Overwhelming Evidence 117.41* [105.48, 130.41] 1.000 >100 Overwhelming Evidence 3.70* [ 3.48, 3.91] 1.000 >100 Overwhelming Evidence NA NA NA NA NA 0.29* [0.16, 0.50] 1.000 >100 Overwhelming Evidence
Hurdle Intercept 1.21 [ 0.87, 1.69] 0.868 0.012 Very Strong Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Conditional Within-Person Effects
Daily persuasion experienced 1.03 [ 0.97, 1.08] 0.828 0.028 Very Strong Evidence for Null 1.03 [ 1.00, 1.06] 0.966 0.036 Strong Evidence for Null 0.00 [-0.04, 0.05] 0.553 0.004 Very Strong Evidence for Null 0.85* [0.71, 0.99] 0.980 0.253 Moderate Evidence for Null 0.84 [0.69, 1.01] 0.966 0.199 Moderate Evidence for Null
Daily persuasion utilized (partner’s view) 1.03 [ 0.98, 1.08] 0.899 0.035 Strong Evidence for Null 1.02 [ 0.99, 1.05] 0.888 0.013 Very Strong Evidence for Null 0.02 [-0.02, 0.07] 0.829 0.007 Very Strong Evidence for Null 1.03 [0.84, 1.24] 0.607 0.041 Strong Evidence for Null 1.13 [0.85, 1.54] 0.789 0.077 Strong Evidence for Null
Daily pressure experienced 0.89* [ 0.80, 0.99] 0.984 0.484 Weak Evidence for Null 0.94 [ 0.88, 1.01] 0.960 0.072 Strong Evidence for Null -0.04 [-0.14, 0.07] 0.767 0.014 Very Strong Evidence for Null 1.84* [1.18, 2.68] 0.994 3.366 Moderate Evidence 2.04* [1.03, 4.57] 0.979 1.451 Weak Evidence
Daily pressure utilized (partner’s view) 0.94 [ 0.86, 1.03] 0.915 0.072 Strong Evidence for Null 0.98 [ 0.92, 1.05] 0.714 0.015 Very Strong Evidence for Null -0.02 [-0.14, 0.08] 0.661 0.011 Very Strong Evidence for Null 1.22 [0.71, 2.06] 0.808 0.149 Moderate Evidence for Null 1.44 [0.58, 4.04] 0.799 0.251 Moderate Evidence for Null
Daily pushing experienced 1.03 [ 0.96, 1.10] 0.775 0.032 Strong Evidence for Null 1.03 [ 0.98, 1.08] 0.900 0.025 Very Strong Evidence for Null 0.02 [-0.04, 0.09] 0.768 0.008 Very Strong Evidence for Null 1.17 [0.97, 1.43] 0.946 0.145 Moderate Evidence for Null 1.28* [1.01, 1.64] 0.979 0.395 Weak Evidence for Null
Daily pushing utilized (partner’s view) 0.99 [ 0.93, 1.05] 0.618 0.021 Very Strong Evidence for Null 1.02 [ 0.97, 1.06] 0.771 0.011 Very Strong Evidence for Null 0.08* [ 0.01, 0.14] 0.984 0.079 Strong Evidence for Null 0.91 [0.71, 1.17] 0.770 0.064 Strong Evidence for Null 0.89 [0.60, 1.31] 0.735 0.096 Strong Evidence for Null
Day 1.01 [ 0.89, 1.13] 0.551 0.042 Strong Evidence for Null 0.97 [ 0.91, 1.04] 0.785 0.020 Very Strong Evidence for Null 0.26* [ 0.15, 0.37] 1.000 50.164 Very Strong Evidence 1.47 [0.75, 2.89] 0.871 0.256 Moderate Evidence for Null 1.64 [0.77, 3.51] 0.900 0.364 Weak Evidence for Null
Daily weartime NA NA NA NA NA 1.00* [ 1.00, 1.00] 1.000 0.070 Strong Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Conditional Between-Person Effects
Mean persuasion experienced 1.01 [ 0.74, 1.38] 0.533 0.102 Moderate Evidence for Null 1.10 [ 0.82, 1.46] 0.748 0.075 Strong Evidence for Null 0.33 [-0.21, 0.90] 0.889 0.118 Moderate Evidence for Null 1.12 [0.40, 3.11] 0.586 0.204 Moderate Evidence for Null 1.99 [0.60, 6.80] 0.870 0.454 Weak Evidence for Null
Mean persuasion utilized (partner’s view) 0.98 [ 0.72, 1.33] 0.543 0.106 Moderate Evidence for Null 0.98 [ 0.73, 1.30] 0.559 0.060 Strong Evidence for Null 0.22 [-0.32, 0.79] 0.790 0.076 Strong Evidence for Null 1.38 [0.45, 4.31] 0.711 0.262 Moderate Evidence for Null 1.90 [0.52, 7.15] 0.829 0.417 Weak Evidence for Null
Mean pressure experienced 1.15 [ 0.80, 1.65] 0.772 0.158 Moderate Evidence for Null 0.98 [ 0.73, 1.31] 0.560 0.062 Strong Evidence for Null -0.31 [-0.86, 0.23] 0.874 0.101 Moderate Evidence for Null 3.51* [1.18, 10.71] 0.990 3.142 Moderate Evidence 17.91* [2.33, 159.59] 0.997 20.152 Strong Evidence
Mean pressure utilized (partner’s view) 0.89 [ 0.62, 1.29] 0.744 0.155 Moderate Evidence for Null 0.97 [ 0.73, 1.28] 0.594 0.060 Strong Evidence for Null -0.30 [-0.84, 0.24] 0.871 0.099 Strong Evidence for Null 1.17 [0.37, 3.66] 0.612 0.239 Moderate Evidence for Null 2.27 [0.24, 19.14] 0.769 0.597 Weak Evidence for Null
Mean pushing experienced 1.32 [ 0.84, 2.07] 0.891 0.322 Weak Evidence for Null 0.97 [ 0.65, 1.45] 0.553 0.080 Strong Evidence for Null 0.22 [-0.54, 0.99] 0.714 0.091 Strong Evidence for Null 1.23 [0.28, 5.62] 0.605 0.308 Weak Evidence for Null 0.83 [0.12, 6.30] 0.580 0.405 Weak Evidence for Null
Mean pushing utilized (partner’s view) 1.40 [ 0.88, 2.21] 0.925 0.454 Weak Evidence for Null 1.25 [ 0.83, 1.85] 0.859 0.152 Moderate Evidence for Null 0.35 [-0.39, 1.12] 0.825 0.121 Moderate Evidence for Null 0.11* [0.02, 0.64] 0.992 7.898 Moderate Evidence 0.08* [0.01, 0.66] 0.990 7.042 Moderate Evidence
Mean weartime NA NA NA NA NA 1.00 [ 1.00, 1.00] 0.909 0.000 Very Strong Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hurdle Within-Person Effects
Hu Daily persuasion experienced 0.65* [ 0.57, 0.73] 1.000 >100 Overwhelming Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily persuasion utilized (partner’s view) 0.75* [ 0.67, 0.84] 1.000 >100 Overwhelming Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pressure experienced 1.23 [ 0.88, 1.71] 0.898 0.152 Moderate Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pressure utilized (partner’s view) 0.66* [ 0.42, 0.95] 0.988 0.890 Weak Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pushing experienced 0.58* [ 0.40, 0.78] 1.000 25.195 Strong Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily pushing utilized (partner’s view) 0.54* [ 0.41, 0.69] 1.000 >100 Overwhelming Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Day 1.09 [ 0.85, 1.42] 0.750 0.067 Strong Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Daily weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hurdle Between-Person Effects
Hu Mean persuasion experienced 0.83 [ 0.36, 1.87] 0.680 0.183 Moderate Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean persuasion utilized (partner’s view) 0.83 [ 0.36, 1.89] 0.671 0.173 Moderate Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pressure experienced 3.31* [ 1.36, 8.35] 0.995 6.142 Moderate Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pressure utilized (partner’s view) 1.79 [ 0.74, 4.48] 0.905 0.428 Weak Evidence for Null NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pushing experienced 0.35 [ 0.11, 1.11] 0.962 1.141 Weak Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean pushing utilized (partner’s view) 0.34 [ 0.11, 1.10] 0.964 1.261 Weak Evidence NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Hu Mean weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.32 [0.24, 0.42] NA NA NA 0.31 [0.24, 0.40] NA NA NA 0.60 [0.47, 0.78] NA NA NA 0.82 [0.47, 1.25] NA NA NA 1.18 [0.75, 1.74] NA NA NA
sd(Hurdle Intercept) 0.90 [0.69, 1.17] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Daily persuasion experienced) 0.12 [0.08, 0.17] NA NA NA 0.05 [0.02, 0.08] NA NA NA 0.04 [0.00, 0.10] NA NA NA 0.18 [0.01, 0.43] NA NA NA 0.21 [0.01, 0.51] NA NA NA
sd(Daily persuasion utilized (partner’s view)) 0.09 [0.05, 0.13] NA NA NA 0.06 [0.03, 0.09] NA NA NA 0.08 [0.01, 0.13] NA NA NA 0.22 [0.01, 0.51] NA NA NA 0.50 [0.12, 0.97] NA NA NA
sd(Daily pressure experienced) 0.08 [0.00, 0.24] NA NA NA 0.05 [0.00, 0.14] NA NA NA 0.08 [0.00, 0.24] NA NA NA 0.56 [0.09, 1.14] NA NA NA 1.12 [0.16, 2.43] NA NA NA
sd(Daily pressure utilized (partner’s view)) 0.07 [0.00, 0.19] NA NA NA 0.04 [0.00, 0.12] NA NA NA 0.09 [0.00, 0.26] NA NA NA 0.50 [0.02, 1.55] NA NA NA 0.97 [0.04, 2.73] NA NA NA
sd(Daily pushing experienced) 0.11 [0.04, 0.19] NA NA NA 0.07 [0.01, 0.15] NA NA NA 0.06 [0.00, 0.14] NA NA NA 0.22 [0.01, 0.50] NA NA NA 0.25 [0.01, 0.60] NA NA NA
sd(Daily pushing utilized (partner’s view)) 0.09 [0.02, 0.17] NA NA NA 0.04 [0.00, 0.10] NA NA NA 0.07 [0.00, 0.17] NA NA NA 0.18 [0.01, 0.57] NA NA NA 0.31 [0.01, 0.95] NA NA NA
sd(Hu Daily persuasion experienced) 0.18 [0.02, 0.34] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily persuasion utilized (partner’s view)) 0.17 [0.03, 0.33] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure experienced) 0.31 [0.01, 0.89] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pressure utilized (partner’s view)) 0.34 [0.01, 0.99] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing experienced) 0.64 [0.32, 1.07] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sd(Hu Daily pushing utilized (partner’s view)) 0.32 [0.05, 0.64] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Additional Parameters
ar[1] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
ma[1] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
cosy NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
nu NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
shape NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sderr NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
sigma 0.68 [0.66, 0.71] NA NA NA 0.57 [0.56, 0.59] NA NA NA 0.96 [0.94, 0.98] NA NA NA NA NA NA NA NA NA NA NA NA NA 
report::report_system()

Analyses were conducted using the R Statistical language (version 4.4.1; R Core Team, 2024) on Windows 11 x64 (build 22635)

report::report_packages()
  • beepr (version 2.0; Bååth R, 2024)
  • R.methodsS3 (version 1.8.2; Bengtsson H, 2003)
  • R.oo (version 1.26.0; Bengtsson H, 2003)
  • R.utils (version 2.12.3; Bengtsson H, 2023)
  • brms (version 2.21.0; Bürkner P, 2017)
  • Rcpp (version 1.0.13; Eddelbuettel D et al., 2024)
  • bayesplot (version 1.11.1; Gabry J, Mahr T, 2024)
  • lubridate (version 1.9.3; Grolemund G, Wickham H, 2011)
  • DHARMa (version 0.4.6; Hartig F, 2022)
  • wbCorr (version 0.1.22; Küng P, 2023)
  • tibble (version 3.2.1; Müller K, Wickham H, 2023)
  • R (version 4.4.1; R Core Team, 2024)
  • openxlsx (version 4.2.7.1; Schauberger P, Walker A, 2024)
  • ggplot2 (version 3.5.1; Wickham H, 2016)
  • forcats (version 1.0.0; Wickham H, 2023)
  • stringr (version 1.5.1; Wickham H, 2023)
  • rvest (version 1.0.4; Wickham H, 2024)
  • tidyverse (version 2.0.0; Wickham H et al., 2019)
  • readxl (version 1.4.3; Wickham H, Bryan J, 2023)
  • dplyr (version 1.1.4; Wickham H et al., 2023)
  • purrr (version 1.0.2; Wickham H, Henry L, 2023)
  • readr (version 2.1.5; Wickham H et al., 2024)
  • xml2 (version 1.3.6; Wickham H et al., 2023)
  • tidyr (version 1.3.1; Wickham H et al., 2024)
  • knitr (version 1.48; Xie Y, 2024)
  • kableExtra (version 1.4.0; Zhu H, 2024)
report::cite_packages()